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ABSTRACT 



■ Here we present a novel N-body simulation technique that allows us to compute en- 
\ semble statistics on a local basis, directly relating halo properties to their environment. 

This is achieved by the use of an ensemble simulation in which the otherwise indepen- 
\ dent realizations share the same fluctuations above a given cut-off scale. This produces 

a constrained ensemble where the LSS is common to all realizations while having an 
Q ■ independent halo population. By generating a large number of semi-independent real- 

I izations we can effectively increase the local halo density by an arbitrary factor thus 

c/) . breaking the fundamental limit of the finite halo density (for a given halo mass range) 

determined by the halo mass function. 

This technique allows us to compute local ensemble statistics of the matter/halo 
distribution at a particular position in space, removing the intrinsic stochasticity in 
the halo formation process and directly relating halo properties to their environment. 
This is a major improvement over global descriptors of the matter/halo distribution 
r ■ which can not resolve local variations. 

OO ' We introduce the Multum In Parvo (MIP) constrained ensemble simulation con- 

sisting of 220 realizations of a 32 h ^ Mpc box with 256^ particles each. This is 

■ equivalent in terms of effective volume and number of particles to a box of ~ 193 

h~^ Mpc of side with ~ 1540^ particles containing ^ 5 x 10^ haloes with a minimum 

CN| . mass of 3.25 x 10^ h~^ M©. The MIP simulation stands apart from all previous N- 

body simulations in its unprecedented high equivalent particle density and local halo 
density. 

We illustrate the potential of the technique presented here by computing the local 
' mass function at several characteristic environments and along a path from the center 

, of a void to its border. We can study for the first time the effect of local environment 

■ in the height, shape and characteristic mass of the halo mass function. 

Future observations will provide detailed models of the galaxy density distribution 

allowing the use of constrained realizations in combination with ensemble simulations 
to derive accurate "mean ensemble" properties of the local Universe. 

Key words: Cosmology: large-scale structure of Universe; methods: data analysis, 
N-body simulations 



1 INTRODUCTION 

N-body simulations are one of our most valuable tools in 
cosmology. They allow us to follow the evolution of cosmic 
structures from the linear regime to the present time over a 
wide range of scales and masses. Current state-of-the-art N- 
body simulations contain billions of particles and can resolve 
from large structures like superclusters and clusters of galax- 
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ies down to subhaloes and dwarf galaxies ([Springel et al.l 
I2OO5I : iDiemand et~aIll2QQ7l : ISpringel et al.ll20Q8l ). 

One area where N-body simulations have been par- 
ticularly useful is in the study of the formation and 
evolution of haloes and their relation with their large- 
scale environment. Haloes are affected by their environ- 
ment through the local densi ty , tida l field, matter ac- 
cretio n , mergers, et c . (I White (Il984ll: Bvrd Valtonenl 
(ll99Qll:[l^ev Cold (ll993ll:IZabludoff Mulch^ (Il998l l: 
Ivan den BoschI (|2QQi l: ICottlober et all (|200lh among many 
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others). In order to properly understand the role of local 
environment in the process of galaxy formation/evolution 
it is crucial to be able to perform local correlations between 
haloes and their local large-scale structure (LSS). An impor- 
tant application of this is the use of simulations to interpret 
real observations of specific (and often complex) environ- 
ments such as our local neighborhood. This is an non-trivial 
task as the many "problems" existing between observations 
and simulations show. Examples of s uch discrepancies in- 
clude the "cold Hubble flow" problem (|Sandage et al.lll97^: 
Sandagd [19861 ). the "local velocity anomaly" (i Tullv etaL 



20081 ) . the exces s of luminous galaxies above the plane of 



the local sheet (^ Peebles Nusserl l2010l) and, to s ome de- 
gree, the "missing satellite" problem ([Klvpin et al 1 IT999 ). 

Until recently our current set of theoretical and numer- 
ical analysis tools was limited on global descriptors such 
as the two-p oint correlation function an d its generalization 
to n— points (IPeeblesI 119801 : ISzapudil [T998). marked correla- 
tions (IShethll2005l\ topology via the Minkowski functionals 
(jMecke et al .111 9941 : ISchmalzing et al]ll999l \ etc. These tools 
provide a good description of the global clustering and scal- 
ing relations of the matter in the Universe. However, since 
such measurements contain no phase information they are 
insensitive to local variations that may arise in the diver- 
sity of environments in the Cosmic Web. Some important 
advances have been done in this respect with the introduc- 
tion of LSS characterization algo rithms based properties o f 



the LSS such as the density field (lAragon-Calvo et al 
and related Hessia n-based de scriptors (IZhang et al 



ICautun et alll2012ll t idal field (iHahn et al.ll2007l). shear ten- 
sor ( Libeski nd et al.1 '2012). topologv (^ Sousbie et all l2008l : 
lArag on-Ca lvo et al.| [2010) among many others. Using such 
tools haloes can be assigned to particular environments and 
then study the dependence of their properties. However, af- 
ter the assignment of haloes of specific environments, invari- 
ably the halo sample is integrated into global environment- 
based samples, thus losing the locality gained by the LSS 
characterization. The main reason for this is the need for a 
large halo sample in order to perform meaningful statistical 
analysis. This would be in principle sufficient if all envi- 
ronments could be completely characterized using the algo- 
rithms mentioned above. However, even the most complete 
characterization tools available can not fully describe the 
complexity of environments found in the Cosmic Web. 



1.1 Finite halo sampling 

The problem of correlating halo properties with their par- 
ticular environment can be then described as a "sampling 
problem" . Local studies of haloes and their environment are 
ultimately limited by the local density of haloes which de- 
pends on the halo mass given by the halo mass function 
(jPress k Sc hechter 1974). The halo mass function sets lim- 
its to the number of haloes of a given mass that we expect 
to find in a given environment. We refer to this as the spe- 
cific halo number density. As the halo mass increases their 
specific number density decreases until there are not enough 
haloes to perform any statistical analysis or even until there 
are no ha loes at all as in the cas e of the regions at the center 
of voids ([Cottlober et al.ll2003l V In the Universe and stan- 
dard computer simulations we only have one realization to 
make measurements. At any point in space and at a given 



halo mass range there is a predefined halo number density 
given by the halo mass function. The specific halo density 
does not depend on the mass resolution of our simulation. A 
higher density of particles will only increase the number of 
lower-mass haloes. This sets fundamental limits to the kind 
of analysis we can perform on a local basis. One could for 
instance measure statistics of haloes inside a sampling win- 
dow small enough to be sensitive to local variations. While 
this may work for very limited cases such as low-mass haloes 
and high density regions, in general one simply does not have 
enough haloes to obtain statistically significant results. This 
effect is particularly severe in low density regions such as 
voids and walls where the density of halos massive enough 
to host luminous galaxies is of the order of one halo per 
several cubic megaparsec. 

In this paper we present a new technique that solves 
the problem of a finite specific halo density. The technique is 
based on a constrained ensemble simulation where all semi- 
independent realizations share the same large-scale fluctu- 
ations. Ensemble simulations are becoming an important 
tool to study the statistical properties of the matter/galaxy 
distribution where a l arg e sampling v olume is re quired 
^ leiksin k Whitelll999l: iTakahashi et al.ll2 009: Sm ithll2009l: 
Schneider et al.1 l201ll: iForero- Romero et al.. „201inOrbanl 
20121 ). For instance. ITakahashi et al.1 (|2009h ran 5000 inde- 



pendent realizations of a 1 h ^ Gpc box with 256^ par- 
ticles ea£h__to_stu^yth^^ Acoustic Oscillations (see 



also I Schneider et al 



20071) 

]2009( : iForero-Romero et al 



At smaller physical scales, 
(|2011r ) used ensemble simulations from 
constrained random realizations of the local supercluster in 
order to identify structures of interest. These authors ran 
200 constrained realizations of the local supercluster and 
identified "local group" systems in only 3 realizations. From 
these 3 realizations higher resolution simulations were gen- 
erated and analyzed. As far as we know the remaining 197 
realizations in the ensemble were not used. 



2 LOCAL ENSEMBLE SIMULATIONS 

The overall structure and dynamics of the Cosmic Web is 
mainly determined by large-scale fluctuations with smaller 
fluctuations playin g only a minor role i n the overall picture 
(iLittle et al.1 Il99ll : ISuhhonenko et al.1 l201ll : lEinasto et aD 
l2011al ). Since the primordial density field is a linear combi- 
nation of fluctuations it can be separated into two regimes: 
one responsible for shaping the large-scale features of the 
Cosmic Web and another with smaller galaxy-size fluctua- 
tions as: 



fl) 



where 5 is the total primordial density fleld, ^(< /Ccut) is the 
large-scale fluctuations fleld: 

5{<K^,) = 5{k)W{k). (2) 
and kcut) is the high-pass flltered density fleld given by: 
k.^t) = S{k) {1-W{k)). (3) 
Where W{k) is a fllter with a cutoff frequency kcut given by: 
0, 



W{k) = 



k ^ kcut 
k ^ kcut 



(4) 
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Following equation [T] we can generate multiple semi- 
independent realizations from the same primordial density 
field by fixing the large-scale fluctuations S{< kcut) be- 
tween realization and allowing the small-scale fluctuations 
^(^ kcut) to vary. The choice of kcut determines the scale 
and equivalent mass below which realizations are indepen- 
dent between each other: 

Density fluctuations below this scale are statistically inde- 
pendent at the initial conditions and latter become corre- 
lated by the power transfer from large scales. These semi- 
independent realizations define multiple paths for initially 
statistically independent small-scale objects that neverthe- 
less converge to the same large-scale patterns by their shared 
large-scale fluctuations. Haloes with masses larger than 
M(kcut) are not independent between relations offering a 
unique opportunity to study their "persistent "properties 
across the ensemble. These properties reflect the large-scale 
environment of the halo as the small-scale contributions are 
averaged-out by the ensemble. 

2.1 Nested ensembles 

The idea of a constrained ensemble can be extended to it- 
eratively generate nested constrained ensembles where for 
each realization i\ 

6' = 6{<kcut) + S{^kcuty (6) 

we generate a new ensemble defined by a new cut-off scale 
A:cut2 where A;cut2 > kcut • Realizations in the nested ensemble 
will share structures at scales larger than /Ccut2 and will be 
independent below those scales: 

S2 = ^(< kcut) + S{'^ kcut, < kcut^y + ^(^ A:cut2). (7) 

Here the subscript denotes the level of recursion in the en- 
semble. 

Designing a nested constrained ensemble one could de- 
fine a cutoff scale kcut corresponding to cluster-group size 
structures and having a galaxy-size independent halo popu- 
lation between realizations. The next recursion level would 
be then determined by a higher A;cut2 corresponding to 
Milky- Way size galaxies thus generating a sub-halo indepen- 
dent population between realizations for every galaxy-size 
object in the first level of the ensemble. 

The number of realizations grows with the number of 
nested ensembles as where N is the number of realiza- 
tions per ensemble and m is the level of recursion. This expo- 
nential behaviour makes it unfeasible to generate more than 
a just few recursion levels. In the present work we present a 
one- level ensemble. A two- level nested ensemble is in prepa- 
ration at the time of writing. 

2.2 Practical advantages 

The self-contained nature of each realization offers several 
practical advantages over standard N-body simulations: 

• Trivial to parallelize (run and post-procesing) . Realiza- 
tions are self-contained and can be run and analised on an 
independent pipeline. 



• Computing time increases linearly with the number of 
realizations. 

• Adding more realizations is trivial as it only involves 
generating new small-scale fluctuations. 

• No need for custom read-write routines as in the case 
of standard massive simulations where parallel 10 routines 
are necessary given the shear size of the datasets. 

• Efficient storage of many relatively small snapshots. 

It may seem that the above advantages are the same as 
if we simply ran a large number of independent realizations. 
However, in our case each realization corresponds to the 
same LSS configuration so all realization can be effectively 
considered to be the same simulation with different halo 
formation paths. 



3 LOCAL ENSEMBLE STATISTICS 

The concept of ensemble implies the existence of statistically 
independent events. The ensemble average of a variable of 
the system is given by 

N 

M ^ Wipi (8) 

i=l 

where Wi is the measured variable at event i and pi is the 
probability of observing the event i with pi — 1/N when 
all events have the same probability of being "observed". 
Angled brackets denote the ensemble average. In general it is 
desired to know the moments of a variable with respect to its 
ensemble average. The n^^ central moments of a distribution 
are then: 

{\w - {w)n = lim -(«,»". (9) 

where (w) is the ensemble average defined in equation [8] 
and N is the number of realizations or observations in the 
ensemble. It is common to invoke ergodicity and interchange 
volume averages by ensemble averages. In our case however, 
since realizations in the ensemble are correlated above kcut 
this is not necessary and we can define the local ensemble 
moments at position x as: 

(|«;(x)-Hx))r>= lim lVWx).-(«;(x)»". (10) 

N^OG iV ^ ^ 

i=l 

In practice we compute the local n^^ central ensemble mo- 
ment by considering all M sampling points (halos or parti- 
cles) inside a sampling window of radius R located at posi- 
tion Xi across all N realizations: 

N M 

{\w{^)-{w{^))nnn = ^ E 7^ ^(^(x.,,)-(^(x))h)^(11) 

where w{:x.ij) is the variable measured at sampling point Xi 
inside the sphere corresponding to realization j and {w{'x.))r 
is the local ensemble mean. 

3.1 Realization stacking and signal-to- noise 

The cross-correlated nature of the realizations in the con- 
strained ensemble allows us to "stack" them into one single 
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LIMITED BY COMPUTATIONAL RESOURCES | 
■ 

AQUARIUS 



k [h Mpc"'] 

Figure 1. Cross-correlation coefficient between the initial condi- 
tions of two realizations in the ensemble. 



field. By doing so we 'average-out" the effect of stochastic 
processes and focus on the effect of the environment defined 
by scales larger than the cut-off scale /Ccut • This is a similar 
process than the "image stacking" technique used in image 
processing for noise reduction. The signal-to-noise ratio of 
the stacked ensemble is given by: 

s,n=^. (12) 

Where N is the number of realizations in the ensemble. Al- 
though in practice we can not create an infinite number of 
realizations one can generate enough realizations in order to 
increase the signal-to-noise ratio to any desired level. 




Halo mass 

Figure 2. Comparison of the number of particles per halo mass 
between the MIP and several simulations in the literature. Lines 
indicate simulations with a "continuous" halo population and iso- 
lated symbols indicate a resimulation of individual halos. The 
single-realization MIP is shown as a black line. The MIP ensem- 
ble simulation is shown in red. The lower regime corresponds ha- 
los with masses below M(/ccut)- In this mass regime haloes are 
independent between realizations and can not be stacked into an 
"ensemble halo" . The upper regime of the MIP simulation corre- 
sponds to haloes with masses larger than M(fccut)- In this mass 
regime haloes can be stacked as they are cross-correalted across 
the ensemble. Note that this diagram does not refer to the specific 
halo density (see text for details) which is 220 times larger for the 
MIP compared to standard simulations. The MPI simulation oc- 
cupies a unique region in the diagram in terms of its two-regime 
behaviour and unprecedented particle count for a continuous halo 
population. It has the highest number of particles per halo above 
M„,t. 



4 THE MIP SIMULATION: SETTING UP 

Here we describe the Multum In Parv^ (MIP) ensemble 
simulation. The MIP simulation is an undergoing project 
consisting (in its first stage) of 512 realizations of a 32 
/i~^Mpc box sampled on a 256^ regular grid giving a mass 
per particle of 1.62 x lO^Moh"^. The cosmology corre- 
sponds to the standard ACDM with parameters = 0.3, 
Qa = 0.7, h=0.73, (78 = 0.84 and spectral index n — 0.93. 

The results presented here correspond to 220 out of a 
total of 512 realizations that have been completely run and 
analyzed to identify FoF haloes in all snapshots and Sub- 
Find subhaloes at the final snapshot as described in Section 
14.21 The full MIP ensemble simulation is equivalent in terms 
of volume and number of particles to a 256 /i~^Mpc box 
standard simulation with 2048^ particles. The current sta- 
tus of 220 realization is equivalent to a 32 x (220)^/^ - 193 
h"^Mpc box with 256 x (220)^/^ - 1540^ particles. With 220 
realizations the signal-to- noise ratio is s/n = V220 15. A 
comparison between the properties of a single realization 
and the ensemble simulation is presented in Table [1] 



^ Latin for "much in little". The expression is also used in com- 
puter graphics to denote image pyramids (mipmaps) for texture 
rendering. 



4.1 Generating initial conditions 

The first step in making the constrained ensemble is the 
creation of a "template" initial conditions field from which 
all realizations will be generated. The template field was 
generated using the publicly available grafic code. From 
the template each of the semi-indenpendent realizations was 
generated as described in section [2] (see also Figure IA1|) . In 
practice we do not generate a new random field for each 
realization but instead shift ^(> /Ccut) in real space along 
each dimension by a distance corresponding to /Ccut • For our 
study we set /Ccut = 4 h~^Mpc corresponding to a volume 
containing a mass of ^ 2x 10^^ h~^ M©. For the 32 h~^ Mpc 
box this gives 32/4 = 8 independent shifts per dimension 
for a total of8x8x8 = 512 realizations. From the set of 
initial condition files we then generated Gadget files starting 
at ^ ^ 80 follo wing; the standard Zel'dovich prescription 
(|Zeldovichlll970l l. 

Figure [T] shows the cross correlation coefficient between 
two initial density field realizations computed as: 



(13) 



There is a perfect cross-correlation between the two realiza- 
tions at scales < /Ccut as expected since by construction both 
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z=5 z=l z=0 z=0 (zoom) 




Figure 3. Ensemble simulation. Four realizations are shown at redshifts z = 5, 1 and (first three panels from left to right). The 
stacking of the 220 realizations in the ensemble is shown in the bottom panels. The colors indicate the local density on a thin slice across 
the simulation box. A zoom region of 8 x 8 h"-*^ Mpc centered on a cluster is shown in the right panels. 
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initial density fields share the same large-scale fluctuations. 
At scales {k > kcut) there is practically no cross-correlation 
although there are small noisy "ripples" possibly arising as a 
result of having the same small scales only shifted instead of 
being completely independent realizations. Figure [1] shows 
that while small-scale fluctuations from adjacent realizations 
produced by shifting are not strictly independent from each 
other they are sufficiently separated to minimize significant 
cross-correlation. 



4- 1-1 Ensemble run 

Each realization was ran on 16 processors on the Homewood 
High Performance Cluster (HHPC) at JHU. We stored 50 
snapshots separated in logarithmic intervals in scale factor 
starting at z = 5. The realizations were ran and processed 
in batches of 20 runs in order to keep a low load on the 
HHPC cluster. All the snapshots from the 220 realizations 
plus halo/subhalo catalogues occupy approximately 5.5 TB. 
Each "ensemble snapshot" occupies ^ 100 GB in contrast 
to the much smaller "realization snapshot" with only ^ 460 
MB. 



4.2 Halo/subhalo catalogues 

We identified Friends of Friends (FoF) haloes for all 50 snap- 
shots of each realization using a linking parameter ofb = 0.2. 
Only haloes with a minimum mass of 20 particles were used. 
For each halo we computed extensive physical properties 
such as mass, radius, inertia tensor, angular momentum, 
'^^max, etc. The details of the computation of halo properties 
will be described in detail in a forthcoming paper. Subhaloes 
were identified onl y at the final snapshot at 2; = with the 
SUBFIND code (Sp ringel et "alir2001) which identifies sub- 
haloes as gravitationally bounded substructures inside FoF 
haloes. We also computed physical properties for the sub- 
haloes and stored both halo and subhalo catalogues in a 
database for efficient retrieval. 

Given the large number of realizations, the creation, 
running and analysis of each of the 512 realizations is con- 
trolled with an automated pipeline. The final result is a set 
of halo/subhalo catalogues and most massive progenitor lists 
for FoF haloes. 



5 THE MIP SIMULATION: FIRST RESULTS 

We now describe the general properties of the MIP simu- 
lation starting with individual realizations, the global mass 
function (per realization and ensemble) and the basic prop- 
erties of the most massive haloes in the simulation. 



4-1.2 MIP compared to standard simulations 

Figure [2] shows a comparison between the MIP and sev- 
eral other simulations in the literature: Carmen (from 
"Las Da mas" , McBrid e et al. 2011, in prep). Mil- 
lennium (Springel et al.' ^5), Millennium II (?) and 
Millennium XXL (Angulo et al. 20lJ). We also show 
three high-resolution zoom resi mulations: Aquarius project 
(ISpringel et al.ll20 08h. Phoenix (|Gao et al.|[2012. ) and Rhap- 
sodv (|Wu et al. 1^012). This is by no means a complete list. 
It is intended to show the current state-of-the-art in N- 
body cosmological simulations. The single-realization MIP 
behaves like a standard simulation. It is, for current stan- 
dards, a small-box medium-resolution simulation. It sits 
between the Millennium and Millennium II simulations in 
terms of particle resolution but it contains a much smaller 
volume. On the other hand, the MIP ensemble simulation 
occupies a unique region in the diagram. It has the highest 
number of particles per halo in its "local ensemble range" de- 
termined by M(kcut) (almost one order of magnitude higher 
than millennium II). Haloes in this regime are persistent 
across the ensemble with small variations between realiza- 
tions allowing us to study their ensemble statistics. Haloes 
less massive than Mcut are independent between realizations 
and can be stacked to increase the specific halo density. In 
our case the density of haloes is 220 times larger than in the 
single-realization case. For comparison, the density of haloes 
in the MIP simulation is almost two times higher than the 
density of particles in the millennium I simulation! Given its 
unique properties, one can not directly compare the MIP en- 
semble simulation to other standard simulations. The MIP 
simulations stands apart from standard simulations in terms 
of particles per halo, specific halo density, statistical prop- 
erties, and computational requirements. 



5.1 Individual realizations 

Figure [3] shows an overall view of the MIP ensemble simu- 
lation. Four individual realizations are shown as well as the 
stacked ensemble. From left to right we show the density 
field at redshifts 2; = 5, 1 and and a zoom into the central 
halo in the slice. All realizations share the same large-scale 
features as expected since they evolved from the same large- 
scale initial fluctuations. Already at 2; = 5 we can identify 
common structures between realizations like the large void 
at the top-left corner and the region that will collapse into 
the largest halo at the center of the slice. The small-scale 
fluctuations on the other hand are different in each realiza- 
tion. This becomes more evident as the simulation evolves 
and the fluctuations collapse into haloes with unique sub- 
structure. At ^ = 1 it is clear that even when the LSS is 
roughly the same between realizations it is composed by a 
unique halo population. This effect can be seen in all envi- 
ronments but it is perhaps more dramatic in walls and voids 
where the halo density is very low and often the same region 
is unevenly sampled between realizations. 

All realizations produce a similar cluster at the cen- 
ter of the slice. However, a detailed view inside the cluster 
shows significant differences at small scales. For instance, 
the central cluster in the first two realizations of Figure [3] 
has two large subhaloes while the last two haloes have three 
large subhaloes. The differences at smaller scales are even 
larger making it impossible to identify common low-mass 
subhaloes between realizations. The halo population corre- 
sponding to these scales is effectively independent between 
realizations. Despite this, low-mass haloes trace the same 
large-scale structures as a result of the transfer of power 
from large scales to small scales via gravitational coupling. 




Figure 4. Thirty realizations of the central cluster in Figure O The SubFind subhaloes are shown as open circles scaled with their 
radius. The background image represents the dark matter density field for the corresponding realization. 



8 Aragon-Calvo M.A. 





Vbox [h-i Mpc]3 


Veq [h-1 Mpc]3 Npart 


Nhalos {Z = 0) 


s/n 


single 




32^ 256^ 


2 X 10^ 


1 


ensemble 


32^ 


193^ 1545^ 


~ 5 X 10^ 


15 



Table 1. MPI ensemble simulation, comparison between single and ensemble (N=200) properties. Vbox corresponds to the volume of the 
simulation box and Veq the equivalent volume of the ensemble. Npart is the total number of dark matter particles. Nhaios is the number 
of haloes above the minimum particle limit of 20 particles corresponding to 3.2 x 10^ h~^ Mq. s/n is the signal-to-noise ratio of the 
stacked ensemble. 




5.2 Stacked ensemble 

The ensemble density field is shown in the bottom panels of 
Figure [S] This field was created by stacking the individual 
realizations in the ensemble. The same idea is used in im- 
age processing in the "image stacking" technique where a 
sequence of images taken from the same subject are added 
in order to reduce the noise level (see Section [3.1|) . In our 
case the stacking of realizations "averages out" random fluc- 
tuations of scales smaller than kcut- 

The stacked density field has some interesting proper- 
ties that make it useful for applications where a smooth 
field is needed. The stacking of realizations acts like a sort 
of low-pass filter producing a smooth density field while at 
the same time retaining anisotropic features, in contrast to 
what we would expect if the field was smoothed with a Gaus- 
sian kernel. Note that in contrast to the single realizations, 
the ensemble density field looks remarkably similar at all 
redshifts. The only changes we can see seem to be due to 
scales larger than kcut as smaller scales are averaged by the 
stacking. The ensemble density field can be used as a robust 
large-scale density field tracer which is not affected by small 
scale fluctuations such as haloes which tend to negatively 
affect structure finder techniques based on the density field. 
A good illustration of this is the thin filament located at the 
bottom-right of the central cluster. In all four realizations 
there are significant differences in the substructures defining 
the filament. The ensemble density field, on the other hand, 
shows a well defined tenuous structure with practically no 
substructure and an almost constant density profile. 

The ensemble density field looks tantalizing similar 



to a Lagrangian-smoothed density field (see ("Little e t alj 
9l|; Melot t Shandarin 1993; Avila- Reese e^al. 2Q0ll : 
Einasto et al.1 l2011bl : lAragon-Calvo et al.1 l20ld l2012l )) 



where the primordial density field is smoothed in the lin- 
ear regime and then evolved. Such Lagrangian-smoothed 
simulations contain only large-scale modes and therefore 
no substructure below the smoothing scales. The top-right 
panel of Figure [8] shows the lagrangian-smoothed version of 
the cluster in Figure (8] Note the similarity with the clus- 
ter in the ensemble density field shown in the same figure 
and in the bottom-right panel of Figure [3] The Lagrangian- 
smoothed density field presents several artifacts from the 
density reconstruction procedure but in general it has the 
same properties than the ensemble density field. However, 
the lagrangian-smoothed version does not contain any small- 
scale structure and therefore can not be directly used to 
study haloes. 



5.3 Ensemble halo population 

Figure |4] shows 30 versions of the central cluster in Figure [3] 
identified in different realizations. The ensemble mean mass 
of the central FoF halo is ^ 10^^ h^^M©, computed as the 
mean mass of the corresponding halo across the ensemble. 
This is the second largest halo in the ensemble. We high- 
light the substructure inside the cluster by the position of 
the subhaloes as open circles scaled with their radius. While 
the megaparsec-scale matter distribution is roughly similar 
in all realizations there is a large variety in the properties 
of the halo/subhalo populations. This can be seen in the 
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Figure 6. Halo mass function. Top: the black solid line cor- 
responds to the halo mass function of the entire ensemble. The 
thin gray lines show the halo mass function computed from 20 
individual realizations. Bottom: ratio between the mass function 
of each realization and the ensemble mass function. 
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Figure 7. Top: Halo mass dispersion for the 15 most massive 
haloes in the simulation. Bottom: relative halo mass dispersion. 



most massive subhalo of the cluster which has a variation in 
its size of almost a factor of two between the most extreme 
versions of the cluster shown in the figure. Each realization 
of the cluster shows a unique formation path with differ- 
ent mass and halo accretion histories, internal dynamics etc. 
However, because in all realizations the same LSS is shared 
one expect the properties of the cluster to be highly corre- 
lated across the ensemble (Araya-Melo & Aragon-Calvo in 
preparation) . 



5.3.1 Ensemble Halo Mass Dispersion 

In order to get a better understanding of the differences be- 
tween halo populations in the ensemble we computed the en- 
semble halo mass dispersion for the 15 most massive haloes 
in the simulation. We do this by first defining a "reference" 
sample of haloes identified from a single realization. Next 
we track the reference haloes across all realizations in the 
ensemble. We track a given halo by placing a search win- 
dow at the halo's position in all realizations and selecting 
the most massive halo inside the search window. The search 
window radius is of the same size as the original halo. Fig- 
ure [7] shows the ensemble mass dispersion for the 15 most 
massive haloes. The relative dispersion is shown in the lower 
panel. We normalized the dispersion by the ensemble mean 
halo mass for each halo. The relative dispersion increases for 
decreasing halo mass. The trend in the relative dispersion 
indicates that haloes less massive than ^ 5 x 10^^ are defined 
by ^(^ kcut) (in good agreement with the mass- volume de- 
fined by /Ccut) and therefore can be considered independent 
between realizations. The relative halo mass dispersion can 
be use as a guide for setting a threshold above which haloes 
can be studied on a case-by-case ( M(/Ccut) ^ 2 x 10^^ h~^ 
M0 in our case). At lower masses the relative mass disper- 
sion is larger than 30%. 

Figure [5] shows a more detailed view of the distribution 
of halo mass across the ensemble for the 5 most massive 
haloes in the simulation. The physical location of the five 
most massive haloes in the simulation box is shown on the 
left panel and their ensemble mass dispersion is shown in 
the right panel. The ensemble mass distribution per halo 
has a roughly log-normal distribution. Note the bump in 
the mass distribution of the second largest halo. This is the 
same cluster depicted in Figures [3] and lU The bump seems 
to be associated to the large subhalo seen in practically all 
realizations coming from at the top-right side of the clus- 
ter as small variations between realizations make the FoF 
algorithm link the large subhalo in some cases and not in 
others. 



5.3.2 Ensemble Halo Mass Function 

Figure [6] shows the mass function computed from 20 indi- 
vidual realizations (gray lines) and the full ensemble (black 
line) . We also show the ratio between the individual realiza- 
tions and the ensemble in the bottom panel. The individual 
mass functions closely follows the ensemble mass function 
with some variations becoming more significant for masses 
larger than 10^^ h~^ M©, at which point the mass function 
is dominated the relatively small volume of the simulation 
box. The "knee" in the mass function at M- 10^^ h"^ M© 
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is a particular characteristic of the original template used to 
create the ensemble. The volume-mass defined by kcut Hmits 
the halo-mass dispersion between realization so we do not 
expect large variations in mass between large haloes. From 
the volume-mass defined by kcut we expect to have roughly 
the same halo population above a mass of - 10^^ h-^ Mpc. 



6 TWO ILLUSTRATIVE ENSEMBLE 
STRUCTURES 

In this section we present two structures that (given their 
scale) are persistent across the ensemble and can be stacked 
into an "ensemble structure". We focus on two cases of en- 
vironment extremes: the second most massive cluster in the 
simulation and a large void next to it. Here we offer a gen- 
eral qualitative description of the halo population in both 
environments. We believe that every case presented here de- 
serves a more detailed study that goes beyond the scope of 
this introductory paper. A more detailed analysis will be 
presented in a forthcoming paper in preparation. 

6.1 Ensemble cluster 

6.1.1 Density field 

Figure [8] shows a comparison between a single realization of 
a cluster (left panel) and the stacked ensemble (right panel). 
On the top panels we show the density field (saturated at ^ > 
100) and on the bottom panels the subhalo population on a 
slice of 2 h~^ Mpc thick. The ensemble field was computed 
as the mean of the density field computed at each voxel 
between all realizations. 

The density field in the single realization presents sev- 
eral noisy substructures traced by subhaloes and small fil- 
aments compared to the more smooth ensemble density 
field. The averaging across the ensemble has a similar result 
as "averaging out" the differences in small scales fluctua- 
tions between realizations and emphasizing only the com- 
mon large-scale structures (see also Figure [3|). For compari- 
son we show the Lagrangian-smoothed density field at a = 2 
h~^ Mpc in the insert at the top left panel. The Lagrangian- 
smoothed field has the same large-scale structures than the 
ensemble average. It even has the large infalling clump at the 
top-right of the cluster. However, the Lagrangian-smoothed 
density field lacks of haloes below the mass equivalent to the 
smoothing scale. The ensemble density field is dominated by 
features arising from scales above kcut but in contrast to the 
Lagrangian-smoothed field it also contains a large popula- 
tion of collapse objects. 

The ensemble density field allows us to remove the ran- 
domness in halo properties arising from stochastic processes 
intrinsic to halo formation (nature) by "averaging out" the 
small-scale contributions and emphasizing the large-scale 
environment (nurture). We can now see what is the "true" 
structure of the cluster. This is the structure defined not 
by the primordial small-scale fluctuations inside the cluster, 
but by the large-scale environment of the cluster. The most 
salient feature of the cluster is the large clump of matter on 
the top-right side. This clump is formed by large subhaloes 
present in practically all realizations. The large subhalo is in- 
falling into the main halo through the filament coming from 



the top-right. It is important to emphasize that the clump 
in the stacked ensemble density field does not correspond to 
an individual halo but to an ensemble of subhaloes. Because 
of that we can not directly compute the clump's properties 
from the densify field as in a standard simulation. Instead 
we would have to identify the subhalo across the ensemble 
and compute its ensemble moments. 

6.1.2 Subhalo population 

The bottom panels of Figure [8] show the subhalo population 
inside the second largest halo in the simulation. We included 
haloes with more than 20 particles corresponding to a mass 
of 3.25 X 10^ h~^ Mq This figure illustrates the advantages 
of using of the ensemble over single simulations. The sub- 
halo population on the single realization barely samples the 
central halo and its surrounding structures. Filaments are 
sampled by haloes but not densely enough to trace the fil- 
ament's shape. We can see some dependence between halo 
size and local density but because of the stochastic nature of 
halo formation we can not draw any significant conclusions 
from the single realizations. On the other hand, the ensem- 
ble simulation shows an astounding level of detail. We can 
now clearly see the distribution of haloes/subhaloes inside 
and around the central halo. Subhaloes in the ensemble fill 
all the space in and around the cluster. Even the tenuous fil- 
aments are densely traced by hundreds of haloes. Consider 
that for every halo in the single realization we expect in 
average 220 haloes in the stacked ensemble. Haloes can be 
seen in all environments, from the dense core of the central 
cluster down to the center of the underdense adjacent voids. 
There are ^ 5 x 10^ haloes in the ensemble giving a mean 
halo density of the order of ^ 150 haloes per cubic mega- 
parsec. This value is even higher in high-density regions such 
as the environment of the cluster in Figure [8] where one can 
easily achieve number densities of a thousand objects per 
cubic megaparsec. 

6.2 Ensemble void 

We now focus on a large void adjacent to the central halo 
in Figure [8l The void can be seen on the top-left side of 
the central cluster the same figure. As in the case of the 
cluster, we describe first the density field and then the halo 
population. 

6.2.1 Density field 

The density field from a single realization and the complete 
ensemble are shown on the top panels of Figure (9] There 
are significant differences between the single realization and 
ensemble density fields. The single realization shows more 
substructure inside the void as a result of its internal small- 
scale fluctuations ^gai. We can see several tenuous filaments 
crossing the void in the single realization density field. The 
stacked ensemble density field on the other hand is very 
smooth as a result of the averaging over the 220 realizations. 
The ensemble averaging allows us to highlight the persis- 
tent properties of the void given by its surrounding LSS and 
not its internal random fluctuations. We can now see vari- 
ations in the density profile of the void. The density profile 
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SINGLE REALIZATION STACKED ENSEMBLE 




Figure 8. Single vs. ensemble simulation of a cluster (left and right panels respectively). Top panels: density across the 
simulation box on a thin slice of 8 h~^ Mpc of side and 2h~^ Mpc of thickness. The density field from a Lagrangian-smoothed realization 
with a = 2 Mpc is shown for comparison in the insert of the top-right panel (see text for details). Bottom: subhaloes identified 
inside the slice. The circles are scaled with the halo radius and colored with the value of the ensemble density field from the top right 
panel. For clarity we only show subhaloes with radius smaller than 500 h"-*^ Kpc. The gray background corresponds to the density field 
for one single realization and the ensemble (left ad right panels respectively). 



in the top-left to bottom-right direction is steeper than in 
the bottom-left to top-right direction. This difference is al- 
ready hinted in the single realization density field but it is 
not possible to determine if it is an intrinsic property of the 
void or the result of its internal random fluctuations. 



6.2.2 Halo distribution 

The distribution of FoF haloes is shown in the bottom panels 
of Figure [H Here we plot the FoF haloes instead of SUB- 
FIND subhaloes as done in the case of the cluster (see Figure 
[8]) in order to illustrate how haloes trace the internal struc- 
ture of the void and not the internal structure of haloes 
themselves. The difference between single and ensemble sim- 
ulations in voids is even more striking than in the case of 
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SINGLE REALIZATION 



STACKED ENSEMBLE 




Figure 9. Single vs. ensemble simulation of a void (left and right panels respectively). Top panels: density across the 
simulation box on a thin slice of 16 h"-*^ Mpc of side and 2 Mpc of thickness. Bottom: FoF haloes identified inside the slice. The 
circles are scaled with the halo radius and colored with the value of the ensemble density field from the top right panel. The gray 
background corresponds to the density field for one single realization and the ensemble (left ad right panels respectively). 



clusters (see Figure [8]). Voids are practically devoid of galax- 
ies at the minimum halo mass of our simulation (3.25 x 10^ 
M0). We can see a tenuous filament going from the top 
of the cluster down to its center where it seems to disappear 
as the filament continues beyond the slice. From the sin- 
gle realization it is not possible to determine if the filament 
is the result of the large-scale environment or smaller local 
fluctuations. The extremely low halo density also prevents 
us from computing local halo densities which would be zero 
for a "fair" density estimation or completely dominated by 



discreteness effects. Given the low or even null halo density 
it is not possible to compute any meaningful halo statistics 
inside the void. This is the typical case of a void with its 
intrinsic "observational" limitations. The finite halo density 
prevents us from doing any type of local analysis. On the 
other hand, the ensemble simulation gives a completely dif- 
ferent view of the same void. Haloes cover the entire volume 
of the void even down at its extremely underdense center. 
The stacking of realizations increases the halo density by a 
factor of 220 and the signal-to-noise increases to s/n = 15. 
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By increasing the number density of haloes are are able to 
compute local densities even in the most underdense envi- 
ronments where one would normally be dominated by noise. 
As we will see in the following section we are able to com- 
pute even more elaborated local statistics such as local mass 
functions even at the center of voids. 



7 MASS DEPENDENCE ON HALO BIASING 

Figure [10] shows the tracing of the Cosmic Web by haloes 
in two mass ranges: 10^° - 10^^ h"^ and 10^^ - 10^^ h"^ 
Mq (left and right panels) for the single and staked en- 
sembles (bottom and middle panels). This illustrates in a 
unique way the biasing of haloes with respect to the un- 
derlying dark matter density field (jKaiserl ll984 V Halo bi- 
asing depends on many factors such as mass, scale an d 
halo propert i es (IMo &; White! (|l996l ): IWechsler et al.1 (l2006l ): 
ISmith et all (|2Q0 7") among others). In general, the limited 
halo density in standard realizations makes it impossible to 
study the biasing at particular environment. The unprece- 
dented specific halo density of the MIP simulation allows us 
to compute matter and halo densities in all environments 
and over a wide halo mass range. 

The difference in the halo distribution between the sin- 
gle and staked ensembles is striking for both halo mass 
ranges. Haloes in the stacked ensemble closely trace the LSS 
in contrast with the single-realization where the traced LSS 
looks shallow and dominated by sampling noise. Halo bi- 
asing has the effect of sharpening the features in the LSS, 
this effect that can be exploited to select LSS structures in 
complex environments (Aragon-Calvo in preparation). This 
sharpening effect can not be solely explained by local density 
as we will see in the following section. 



7.1 Bias-induced LSS sharpening 

Figure [11] shows a dramatic example of the bias-induced 
sharpening of LSS structures. The figure highlights a fil- 
ament extending between two cl usters. This is the same 
kind of strai g ht fila ments found bv lColberg et al.1 (|2005l ) and 
iDolag et al.l (|2006l V The high halo density combined with 
their biased distribution sharpens the boundaries of the in- 
ternal structure of the filament unveiling its complex struc- 
ture. This filament is in fact formed by three different fila- 
ments collapsing along their mutual axis. The face-on view 
of the filament in the right panel of Figure [11] clearly shows 
three individual filaments with similar thickness. This kind 
of structures are only observed in stacked ensembles where 
the small-scale fluctuations are averaged-out and the high 
halo density allows for the local effects of halo biasing to be 
observed. 

The distribution of high-mass haloes (top-right panel) 
is an excellent illustration of the effect of environment on 
local biasing. The large void at the center of the slice does 
not contain any halo even after 220 realizations. Note the 
sharp boundaries of the void, a small change in the local 
density translates into a sharp fall in the halo density pro- 
file. Compare the matter (gray line) and staked ensemble 
halo (black solid line) density profiles shown in the bottom 
panels. This local biasing effect can not be observed in the 



single-realization case where there are not enough haloes in 
this mass range to trace the LSS. 

The density profiles shown in the bottom panel demon- 
strate the superior halo sampling of the stacked ensemble 
compared to the single realization. The density field was re- 
constr ucted from the halo distribution using the DTFE tech- 
nique (|SchaaDll2007l ). Voids traced by haloes in the single 
realization are shallower than in the stacked ensemble where 
the cosmic structures have sharp boundaries. The low den- 
sity of haloes in the single-realization makes it impossible to 
draw any conclusions on the behaviour of the local biasing. 



8 LOCAL HALO MASS FUNCTION 

The unprecedented local halo density of the MIP simulation 
allows us to perform local analysis that are in general com- 
puted taking advantage of ergodicity, exchanging volume av- 
erages over ensemble averages. In this section we present a 
direct application of local ensemble statistics by computing 
the local halo mass function. This is fundamentally different 
from previous studies based on global or environment-based 
halo samples. The analysis presented here is performed on 
a point by point basis. 

Figure [12] shows a slice of the density field (left panel) 
and four spherical windows sampling the four basic LSS en- 
vironments: void, wall, filament and cluster. We also sample 
a linear path at 13 positions starting at the center of a void 
and ending at the center of a filament in the boundary of 
the void. The halo mass function was computed from all the 
haloes in the ensemble located inside each sampling window. 
The mass function inside each sampling window is shown 
in the bottom-left panels of Figure 1121 The thick gray line 
corresponds to the ensemble global mass function (see also 
Figure [6]). Sampling inside a void, wall, filament and cluster 
are indicated in both panels by solid, dashed, dotted and 
dashed-dotted lines respectively. 

There is a clear trend in the height, slope and character- 
istic mass of the mass function with LSS environment. Both 
voids and walls have mass functions lower than the global 
mass function. The mass function in filaments and clusters is 
higher than the global mass function. The slope of the mass 
function is also higher for voids co mpared to other environ- 
ments (see Gottlober et al. (2003)). We find that the mass 
function inside the selected filament is basically an enhanced 
version of the global mass function. The mass function at the 
center of the cluster is also higher and its slope is also steeper 
than the global mass function. Around a mass of ^ 10^^ h~^ 
M0 the mass function in the cluster and the filament cross 
and after that point there are more massive haloes in the fil- 
ament than in the cluster. This lack of medium mass haloes 
may the be result of the haloes being inside the cluster and 
therefore not accounted for in the FoF catalogue. A more de- 
tailed study will be presented in a forthcoming paper. The 
mass function inside the wall is a fraction of the general mass 
function but it rapidly drops after 3 x 10^^ h~^ M©. For 
this particular wall the probability of hosting a Milky- Way 
halo is lower than 10~^. Even after 220 realizations there 
are no haloes in that mass range at this particular position. 

The lower-right panel of Figure [12] shows the local en- 
semble mass function sampled along a linear path starting 
at the center of a void and ending at a small filament. The 
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Figure 10. Comparison of the spatial distribution of haloes as function of their mass between a single realization (middle panels) and 
the stacked ensemble (top panels). Two mass ranges are shown in the left and right panels. Haloes are indicated by circles scaled with 
their radius and colored according to their local ensemble density. A detailed view of the halo distribution is shown in a zoom region in 
the top-left panel. The reconstructed density field from the halo distribution is shown in the background with green-red color table. Note 
that there are two color tables shown in the figure. One for the haloes and another for the reconstructed density field. The bottom panels 
show the density profile sampled across the center of the slice indicated by the dashed line on bottom-left panel. Gray lines correspond 
to the ensemble dark matter density profile. Solid and dashed lines correspond to single and ensemble halo density profiles respectively. 
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Figure 11. A filament composed by three merging filaments. The left panel shows a 3D view of a large cluster (red) and a filament next 
to it (blue). The spheres correspond to subhaloes in the stacked ensemble scaled with their radius and colored according to their local 
dynamics (Aragon-Calvo in preparation). The three filaments are indicated by numbers. The right panel shows a "face-on" view of the 
filaments in the volume highlighted in the left panel. We draw circles to indicate the location of the three filaments forming a "cosmic 
fork" . 



local density increases away from the center of the void. The 
value of the smoothed {cfq = 2 h~^ Mpc) density fluctua- 
tions at the two end ends of the path are 82 Mpc = —0.51 
and 82 Mpc = 0.31 respectively. The corresponding raw den- 
sity fluctuation values are 8 — —0.82 and 8 — 2.12. This is 
reflected in the local mass function which also increases its 
height and decreases its slope as we move to denser regions. 
Note how the mass function changes uniformly across the 
sampling points in the void in contrast to the signiflcant dif- 
ferences between environments. The mass function smoothly 
increases indicating a smooth change in the environment. 

In order to quantitatively compare the local mass func- 
tion at several sampling positions we parametrized the halo 
mass function using a Schechter function: 

where 0* determines the height of the curve, a is the slope 
of the power-law scaling in the first part of the curve and 
M* is the characteristic mass deflning the transition point 
to an exponential behavior. 

The mapping between the linear path shown in the top- 
left panel of Figure [12] and the M* — a parameter space is 
shown in the top-right panel of Figure [12] As the sampling 
window moves from the center of the void up to the dense 
fllament the characteristic mass increases and the absolute 
magnitude of the slope a decreases. Note that the linear 
path in configuration space maps into a relatively smooth 
path in the M* — a parameter space. This indicates that in 
low density regions not only one expect less haloes but also 
the ratio between more and less massive haloes decreases. 
A simple scaling of the mass function with local density 
can not explain the change in the slope and characteristic 



mass of the mass function. The increase in the slope of the 
mass function is a strong indication that halo formation is 
suppressed inside voids due to other processes in addition to 
local density. 

We would like to emphasize that the results presented 
here correspond to the local ensemble statistics at the par- 
ticular sampling positions, not just global trends. We can see 
continuous changes in the mass function as we move through 
the simulation box. This is a signiflcant advance with respect 
to previous studies. Here we can, for the flrst time, establish 
a direct link between the halo population and the particular 
properties of the LSS environment at any given position. 



9 CONCLUSIONS AND FUTURE PROSPECTS 

We introduced a novel application of constrained ensem- 
bles where all realizations share the same large-scale fluc- 
tuations while having independent small-scale fluctuations. 
This technique represents a signiflcant improvement over 
standard single-realization N-body simulations. By creat- 
ing a large set of semi-independent realizations and stacking 
them we break the fundamental limit of speciflc halo density 
and are able to trace all LSS environments in great detail, 
from the dense clusters down to the underdense center of 
voids. 

We presented the Multim In Parvo (MIP) simula- 
tion/project consisting on a constrained ensemble simula- 
tion with 512 realizations from which 220 have been fully 
ran and analyzed. The MIP simulation in its current sta- 
tus is equivalent in terms of volume and number of particles 
to a 193 h~^ Mpc box containing ^ 1540^ particles and 
5 X 10^ haloes with a minimum halo mass of 3.25 x 10^ 
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Figure 12. Local ensemble halo mass function The left panel shows a slice of the density field and several sampling windows 
corresponding to the four basic LSS environments: cluster, filament, wall and void (solid, dashed, dotted and dotted dashed lines 
respectively) and 13 samplings going from the center of a void to a filament at its boundary (purple to red color lines). The bottom-left 
panel shows the corresponding local halo mass function computed from the halo ensemble at the four environments. The bottom-right 
panel shows the local mass function computed along the 13 sampling points shown in the top-left panel. We indicate the density field 
smoothed at cr = 2 h"-*^ Mpc at the center of each sampling window. The top-right panel shows the effect of the local environment on 
the mass function parametrized by its slope a and characteristic mass M* (see text for details). The colored filled circles joined by lines 
correspond to the sampling path shown in the top-left panel. 



Mq. As a comparison, the unprecedented halo density 
in the MIP simulation is larger than the particle density of 
the millennium simulation by almost a factor of 2! 

We showed preliminar results on the mass dependence 
on the halo biasing and how this effect can be used to iden- 
tify structures in complex environments. 

The unprecedented specific halo density of the MIP sim- 
ulation allows us to compute ensemble statistics on a local 
basis. We presented the first examples of the local mass func- 
tion computed at several representative LSS environments. 
We also computed the local halo mass function along a path 
through the simulation box and showed the mapping of the 



path in configuration space to a path in the M* — a param- 
eter space. 



9.1 Future prospects 

Given its broad scope this technique can be applied to a 
variety of physical phenomena over a wide range of scales. 
One promising application is the use of stacked constrained 
ensembles of the local supercluster to compute local ensem- 
ble statistics of our own cosmic neighborhood and use this 
information to interpret local observations. 
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In a following series of papers we will explore the effect 
of the LSS on halo properties and the ensemble properties 
of large groups and clusters as function of their LSS envi- 
ronment. 
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Figure Al. Block diagram of the construction of the ensemble 
simulation. 



APPENDIX A: ENSEMBLE GENERATION 

Figure lAll shows the steps involved in the creation of the 
ensemble simulation. We start by defining the cosmological 
model, box size and particle resolution. From that informa- 
tion We then create the template initial conditions file (IC). 
The field is then separated into large and small-scale fluc- 
tuations using the sharp-/c filter defined in equation (4] The 
small-scale fluctuations are the shifted along each dimension 
in configuration space by a physical scale corresponding to 
/Ccut = 4 h~-^ Mpc, producing a total of 8 x 8 x 8 = 512 
independent realizations (another possibility is to generate 
completely independent random small scale fluctuations per 
realization). The large-scale fluctuations are simply repli- 
cated for each realization in the ensemble. Each initial con- 
ditions file is the evolved using the Zel'dovich formalism and 
run. 



